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Abstract We consider a reformulation of QED in which covariant Green functions are used to solve for the 
electromagnetic field in terms of the fermion fields. The resulting modified Hamiltonian contains the photon 
propagator directly. A simple Fock-state variational trial function is used to derive relativistic two-fermion 
equations variationally from the expectation value of the Hamiltonian of the field theory. The interaction 
kernel of the equation is shown to be, in essence, the invariant M. matrix in lowest order. Solutions of the 
two-body equations are presented for muonium like system for small coupling strengths. The results compare 
well with the observed muonium spectrum, as well as that for hydrogen and muonic hydrogen. Anomalous 
magnetic moment effects are discussed. 

1. Introduction 

It has been pointed out in previous publications that various models in Quantum Field Theory (QFT), 
including QED, can be reformulated, using mediating-field Green functions, into a form that is particularly 
convenient for variational calculation [1,2]. This approach was applied to the study of relativistic two-body 
eigenstates in the scalar Yukawa (Wick-Cutkosky) theory [3,4,5]. We shall implement such an approach to 
two-fermion states in QED in this paper. 

The Lagrangian of two fermion fields interacting electromagnetically is (H = c = 1) 

C = ip{x) (i'fdtt - mi - q^A^x)) ip(x) + 4>{x) (i^df, - m 2 - qzj^A^x)) <j>{x) 

-i (d a A p (x) - d A a (x)) (d a A^(x) - d A a (x)) . (1) 
The corresponding Euler-Lagrange equations of motion are the coupled Dirac-Maxwell equations, 

{vfdp - mi) il>(x) = q^A^xWx), (2) 



( l7 ^ - ma) 4>{x) = q 2 ^A^x)^(x), (3) 

and 

d^A v (x) - d v d^(x) = f{x), (4) 

where 

fix) = qi^{x)-f v 'ip(x) + q 2 ~4>{x)Y(t>{x). (5) 

Equations (2)-(3) can be decoupled in part by using the well-known formal solution [6,7] of the Maxwell 
equation (4), namely 

A^x) = A°(x) + J dViVfr " x')fix'), (6) 
where D flv {x — x') is a Green function (or photon propagator in QFT terminology), defined by 

d a d a D^(x - x') - d^d a D av {x - x') = g^8\x - x') 7 (7) 
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and A^(x) is a solution of the homogeneous (or "free field") equation (4) with j^(x) = 0. 

We recall that equation (7) does not define the covariant Green function D^ u {x — x') uniquely. For one 
thing, one can always add a solution of the homogeneous equation (eq. (7) with g^ v — > 0). This allows for 
a certain freedom in the choice of D^ u , as is discussed in standard texts (e.g. ref. [6,7]). In practice, the 
solution of eq. (7), like that of eq. (4), requires a choice of gauge. However, we do not need to specify one 
at this stage. 

Substitution of the formal solution (6) into equations (2) and (3) yields the "partly reduced" equations, 



) V>(*), (8) 



(iVfy - mi) 4>{x) = g l7 " (Al(x) + J d 4 x'D^(x - x')f(x'] 
(i-fdp - ma) c/>(x) = c l2 Y (aI{x) + J d 4 x'D^(x - x')f(x')^ 0(x). 



(9) 



These are coupled nonlinear Dirac equations. To our knowledge no exact (analytic or numeric) solution of 
equations (8)- (9) for classical fields have been reported in the literature, though approximate (perturbative) 
solutions have been discussed by various authors, particularly Barut and his co-workers (see ref. [8,9] and 
citations there in). However, our interest here is in the quantized field theory. 

The partially reduced equations (8)- (9) are derivable from the stationary action principle 

SS [i/j] = 6 J d 4 xC R = (10) 

with the Lagrangian density 

C R = ip(x) (i^dn - nix- qu^A^x)) ip(x) + ~<j>(x) (i^d^ - m 2 - q 2 j f ,A£(x)) <p{x) 

-\ J d±x'f{x')D, u {x-x')f{x) (11) 

provided that the Green function is symmetric in the sense that 

Dn„(x - x 1 ) = D^(x' - x) and D^(x - x') = D vll {x - x 1 ). (12) 

The interaction part of (11) has a somewhat modified structure from that of the usual formulation of QED. 
Thus, there are two interaction terms. The last term of (11) is a "current-current" interaction which contains 
the photon propagator sandwiched between the fcrmionic currents. We shall use this modified formulation 
together with a variational approach to obtain relativistic few-fermion equations, and to study their bound 
state solutions. 



2. Hamiltonian 

We shall consider the quantized theory in the equal-time formalism. To this end we write down the 
Hamiltonian density corresponding to the Lagrangian (11), namely 



where 



Tin = Ho+Ti-i + Ti-ii, (13) 



Ho = ^(x)(-i~a -V + m 1 0)ip(x) +<p t (x) (-i~a ■ V + m 2 (3) <f>(x), (14) 
H T = \ J d 4 x'j»(x')D^(x - x')f(x), (15) 
Tin = q 1 '^(x)-f^A^(x)ip(x) +q 2 ^{x)"f ll A^(x)(l)(x), (16) 
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and where we have suppressed the kinetic-energy term of the free photon field. We construct a quantized 
theory by imposing equal-time anticommutation rules for the fermion fields, namely 

{V Q (x,i),^(y,i)} = {0 a (x,t),^(y,t)} = 5 afj 5 3 (x - y) , 

and all other vanish. In addition, there are the usual commutation rules for the Aq field, and commutation 
of the Aq field operators with the ip and </> field operators. 

To specify our notation, we quote the Fourier decomposition of the field operators, namely 

g J (2ir) \ w p / 
with p = p m = (w p , p), cj p = ^/mf + p 2 and 

1 /2 

= E / (5) t SpsC/ (p ' s) e ~ w ' x + D * sV (p ' s) eW ' x ^ ' (18) 

with p = p m = (^p, p), ^ P = \/ m 2 + P 2 - Note that the mass-mi free-particle Dirac spinors u (p, s), w (p, s) 
where (j^p^ — mi) u (p, s) — 0, (7^p M + mi) w (p, s) = 0, are normalized such that 

u ] (p, s) u (p, cr) = u f (p, s) v (p, ct) = — £g CT , (19) 

?7li 

u 1 * (p, s) u (p, cr) = w f (p, s) u (p, cr) = 0. (20) 

Analogous properties apply to the mass-ra2 spinors U, V. The creation and destruction operators b\ b of the 
(free) particles of mass mi, and e?t, d for the corresponding antiparticles, satisfy the usual anticommutation 
relations. The non-vanishing ones are 

{6 ps , fet CT } = { dps , dt ff } = 5 SCT 5 3 (p - q) . (21) 

Again, the analogous properties apply to the mass-W2 operators B\ B, D\ D. As a concrete example, 
we can think of the mass-mi particles as electrons, and the mass-7712 particles as muons, through any pairs 
of charged fermions could be considered. 



3. Variational principle, two-fermion trial states and equations 

The Hamiltonian formalism of QFT is based on the covariant eigenvalue equation 

P" |V) = Q 1 * |V> , (22) 

where P^ — (^H,P^j is the energy momentum operator, = (E,P) is the 4- vector of energy and momen- 
tum, with E 2 — P 2 =M 2 , where M is the invariant mass of the system. There are very few problem in QFT 
for which exact solution of the (j, = equation (22) can be written down. In practice, it is necessary to use 
approximation methods to solve it, such as the widely used covariant perturbation theory or lattice methods. 
Here we concentrate on the variational approach, which is based on the variational principle 

5 (VI H - E k/>) t=0 = 0, (23) 

which we shall consider in the rest frame i.e. P =0. Variational solutions are only as good as the trial states 
that are used. Thus, it is important that the trial states possess as many features of the exact solution as 



3 



possible. For a system like n + e , the simplest Fock-space trial state that can be written down in the rest 
frame is 

Atrial) = £ / d 3 pF SlS2 (pp psi Dl pS2 |0) , (24) 

S1S2 J 

where F SlS2 are four adjustable functions. We use this trial state to evaluate the matrix elements needed to 
implement the variational principle (23), that is 

(Vwl : H - E : \^ trial ) = £ f <fpF* slS2 (p)^ Sl , 2 (p) (u p + Q P - E) (25) 

S1S2 J 

and 

// \ 5 \i \ QiQ2'm 1 m2 ^ /* d 3 pd 3 q 

{i>triai\-Hi:\4> trlal ) = „ 2^ / / n =^1^2 (P)^i^(q) x 

xu(p, si) 7 "u (q, <Ti) D„„(p - g)7 (-q, <7 2 ) 7^ (-p,* 2 ) , (26) 
where the Fourier transform of the Green function was used: 

D,Ax x') = J J0^D,Ak)e-^-*'\ (27) 

For a particle-antiparticle system like positronium an additional virtual-annihilation interaction term, 

u(p, si) 7^ (-p,s 2 ) (ujp)v (-q, cr 2 ) 7^ (q, ct i) (28) 

appears in (26) [10]. Note that we have normal-order the entire Hamiltonian, since this circumvents the 
need for mass renormalization which would otherwise arise. Not that there is a difficulty with handling mass 
renormalization in the present formalism (as shown in various earlier papers; see, for example, ref. [11] and 
citations therein). It is simply that we are not interested in mass renormalization here, since it has no effect 
on the two-body bound state energies that we obtain in this paper. Furthermore, the approximate trial 
state (24), which we use in this work, is incapable of sampling loop effects. Thus, the normal-ordering of 
the entire Hamiltonian does not "sweep under the carpet" loop renormalization effects, since none arise at 
the present level of approximation. Note, also, that (tp tria i\ ■ Hjj : \ip tria i) = 0, that is the variational trial 
state (24) is insensitive to that part of the interaction Hamiltonian which is linear in Aq (x). This means 
that, with the simple anzatz (24) only stable bound states and elastic scattering can be described, but not 
processes that involve radiation. 

The variational principle (23) leads to the following equation 



]T f d 3 P k + n p - e) f SiS2 ( p )sf: iS2 ( P ) 

S1S2 J 



n, s ,n, ^ [^ F ^ H - i) MZ 2aia2 (p,< l )SF; iS2 (p)=0, (29) 



(2tt) 3 „£i , J Up^q 



where M.° T ^ 2l y ia2 (p,q) is the usual invariant matrix element corresponding to the one-photon exchange 
Feynman diagram: 

M7 t t 2aia2 (p, q) = -u (p, Sl ) (-iqan u (q, <n) *£>„„(? - q)V (-q, a 2 ) (-i<&7") V (-p, s 2 ) , (30) 

As mentioned above, for a fermion-antifcrmion system like positronium we obtain [10] the additional virtual- 
annihilation term (qi = q 2 = e) 

M<$™„^ 2 (p, q) - u (p, Sl ) (-ierf) v (-p, s 2 ) iD^ (w p ) v (-q, a 2 ) {-ierf) u (q, a,) . (31) 
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Note that the Al-matrix arises naturally in this formalism, that is M is not put in by hand, nor does its 
derivation require additional Fock-space terms in the variational trial state (24), as is the case in traditional 
formulations (e.g. [12], [13]). 

In the non-relativistic limit, the functions F S1S2 can be written as 

F SlS2 (p) =J F(p)A SlS2 , (32) 

where the non-zero elements of Ay for total spin singlet (S = 0) states are A 12 = — A 2 i = while for the 
spin triplet (S = 1) states the non-zero elements are An = 1 for m s = +1, A12 = A21 = for m s = 0, 
and A 2 2 = 1 for m s = —1. We use the notation that the subscripts 1 and 2 of A correspond to m s = 1/2 
and m s = —1/2 (or f and j) respectively. Substituting (32) into (29), multiplying the result by A SlS2 and 
summing over s\ and S2, gives the equation 

Hp) + fi(p) - E)F(p) = J d 3 q/C(p, q)F(q), (33) 

where 

£(p,q) = -« E A siS2-M SlS2(Tia2 (p,ci)A aia2 . (34) 

s 1 s 2 a- 1 (T2 

To lowest-order in (|p| , |q|) / (mi,m 2 ) (i.e. in the non-relativistic limit), the kernel (34) reduces to K = 
Q1Q2/ |p — q| 2 , and so (33) reduces to the (momentum-space) Schrodinger equation 

(sH^-^/'vV* - (35) 

where e = E~M and m r = mim 2 /M, M = mi + m2. This verifies that the relativistic two-fermion equation 
(29) has the required non-relativistic limit. 



4. Partial-wave decomposition and classification of states 



In the relativistic case we shall not complete the variational procedure in (29) at this stage to get final 
equations for the four functions F SlS21 because they are not independent in general. We require that the 
trial state must be an eigenstate of the relativistic total angular momentum operator, its projection, and 
parity. Namely 



" J 2 " 




' J(J + 1) " 


J 3 


H 1 trial) " 


mj 


V 




P 



Wtrial) ■ 



(36) 



For a system like positronium charge conjugation invariance is an additional requirement, that is 

C\e+e-)=C\e+e-), (37) 

however this does not apply for mi 7^ m 2 . Explicit forms for the operators J 2 , J3 are given in Appendix A. 
The functions F slS2 (p) can be written in the general form 



£ S1S2 ra slS2 

(p) are the usual spherical harmonics. Here and henceforth we will use the notation p = |p| 
etc. (four-vectors will be written as p M ). The orbital indices ^ SlS2 and m SlS2 depend on the spin indices 
Si and s 2 and are specified by equations (36). The radial coefficients of expansion (38) /; 1 1 2 2 132 (P) also 



where Y p S±S2 
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depend on the spin variables. Substitution of (38) into (24) and then into (36) leads to two categories of 
relations among the adjustable functions. 

Mixed-spin states 

In this case £ SlS2 = I = J and the general solution of the system (36) is 

F SlS2 (p) = d4'»)(p) + C 2 i^(p), (39) 

where C\ and C 2 are arbitrary constants. F^sl(p) and Falsi (p) are functions, which correspond to pure 
singlet states with the total spin S = and triplet states with S = 1 respectively. The singlet functions have 
the form 

4 S f 2 ) (P) = /i 1 t'(P)^ mS152 (P) ) (40) 

where mn = m 2 2 = and m\i = m 2 i = mj. The relations between f[^ J {p) and f^^ J (p) involve the 
Clebsch-Gordon (C-G) coefficients C( s fl)™«i«2, that is 

&{p)=C^ m >^f J {p), (41) 

as it shown in Appendix A. We see that the spin and radial variables separate in the sense that the factors 
fsis 2 (p) have a common radial function f J (p), that is for the singlet functions we obtain 

^(P) = C^ m ^f J (p)Y™^ (p). (42) 

The C-G coefficients C^ m '^ have a simple form: C^ mi1 = C {sg ) m22 = 0, C^ sg ) mi2 = -C^ sg ) m21 = 1 
(see Appendix A). Thus the nonzero components of F^"%l(p) arc F^ 9 \p) = F^ 9 \p), F^ 9 \p) = F 2 ^ 9 ^(p). 
The triplet functions have the form 

F^l(p) = f^ 2 J (p)Yj lsiS2 (p), (43) 

where 

mn =mj — 1, mi2 = m 2 i = to,/, to 22 = to,/ + 1. (44) 
The expressions for f£ lS2 (p) involve the C-G coefficients Cj^ 7 ™ 5 for S = 1 listed in Appendix A, that is 

ftsY(p)-C^f( P ), (45) 



where the index m s is defined as 



m s = +1, when m SlS2 — mn, 

m s — 0, when m SlS2 — m i2 = m 2i , (46) 
m s = —1, when m SlS2 — m 22 . 



Thereupon, the triplet function is 

*&(p) = c%; m 'f j (p)Yr-^(p). (47) 

We need to note that (42) is true for J > 0, while (47) is true for J > 1. Thus the coefficient C 2 in (39) is zero 
when J = 0. In other words, for J = 0, only the pure singlet state arises. For a system like positronium the 
requirement (37) decouples the singlet and triplet states for all J. Indeed, the charge conjugation cigenstates 
are 

\ag) = C (ss)m — J d 3 pf J (p)Yj m '(p)bl Si Dl pS2 | 0) (48) 

s l s 2 

with C = (— 1) J for the pure singlet states, and 

M = E C T2 J , msiS2 J d 3 pf J (p)Yr(p)bl Sl Dl pS2 I 0) (49) 
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with C= (-1) J+1 for the pure triplet states, as it discussed in Appendix A. 

The states (48) and (49) diagonalize the Hamiltonian (13) only for particle-antiparticle systems when 
mi = rri2. Thus, for positronium-like systems, the states can be characterized by the spin quantum number 
S, and the mixed states (39) separate into singlet states (parastates S = 0) and triplet states (orthostates 
S = 1). For distinct particles (mi ^ m 2 ) C is not conserved and there is no separation into pure singlet and 
triplet states in general. Thus for arbitrary mass ratios we need to diagonalize the expectation value of the 
Hamiltonian (13) (Appendix E). This can be achieved by the following linear transformation 



where U is the unimodular matrix 



with components 



\sg q ) 

\tr g ) 



U 



U 



\sg) 

\tr) 



a —b 
b a 



where 



£ = (4 ((mi - ma) / (mi + m 2 )f J (J + 1) + l) 



'1/2 



(50) 

(51) 

(52) 
(53) 



The new states, which diagonalize the expectation value of H, shall be called quasi-singlet | sg q ) and quasi- 
triplct I tr q ) states 



\*9q) 



E C ^/ roS1S2 / d 3 pf J (p)Y^(p)bl 1 Dl pS2 \0) 

s l s 2 

tr g ) = E^™/ rosi82 / d 3 pf J (p)Yr(p)blD[ 



(54) 



-PS2 



0), 



where the coefficients (j^f™"!"* an( j Jm * 1S2 are listed in Table 1, and satisfy the following condition 



Jmj 



Jmj 



E ( 



C 



Jmj 



X ( 



c 



(st)Jm v 
Jmj 



2(2J+1). 



(55) 



Table 1. The quasi-state coefficients (j^ Jms i s 2 f QJ . s _ s ^ 



inn = rn.J — 1 



T7T 



mi2 = mj 



m 2 i = mj 



m 2 2 = mj + 1 



, 1 (J+mj)( J-mj + l) 
1 J(J+1) 



T7T 



— + a 



— 



1 / {J-mj)(J+rnj + l) 
~° 1 J(J+1) 



TJT 



s t 



( (J+mj)(J-mj + l) \ 

a \ jv+i) ) 



— + b 



-b 



(,J-rnj)(J+mj 
J(J+1) 



+ 3 



Note that these coefficients differ from C-G coefficients, because of the coupled system we are dealing 
with. Remember that these coupled quasi-states arise only for J > 0. For J = purely 5 = states occur. 
Quasi-singlet and quasi-triplet states are both characterized by the same quantum numbers J, to,/ and 
P = {— 1) J+1 . Because of the unimodularity of matrix (51) we can identify quasi-singlet and quasi-triplet 
states by quasi-spin (like isospin) t = 1/2 with t% = ^1/2, which is a new quantum number (or label). 
However, for our purpose it is more convenient to use the value s = £3 + 1/2, which gives s = s s = or 
s = St = 1 for quasi-singlet and quasi-triplet states respectively. In this case the labels s s and s t reflect 
better the meaning of the indicated quasi-states. It is easy to see from (54) and Table 1 that for positronium 
the quasi states become true singlet (b = 0) and triplet (a = 0) states with different charge conjugation 
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quantum numbers. It is useful to note for subsequent calculations that the coefficients C\ and C2 in (39) are 
C\ = a, C2 = — b for quasi-singlet states s = s s = 0, and C\ — b, C2 = a for quasi-triplet states s — s t = 1. 

The triplet ^-mixing states 

These states occur for l SlS2 = £ = J =F 1 (sec Appendix A). The adjustable functions have the form 

f s1S2 (p) - c^ J - l)m ° f J -\p)Y?:r (p) + c%>y +i)m -f J+ \v)Y?;r m 

where m SlS2 are defined in (46), while the coefficients Cj^^" 7 ^" 1 *, which are precisely the C-G coefficients, 
can be found in Appendix A. Expression (56) involves two radial functions ,f J ^ 1 (p) and f J+1 (p) which 
correspond to I = J — 1 and I = J + 1. This reflects the fact that the orbital angular momentum is not 
conserved and t is not a good quantum number. The system in these states is characterized by J, raj, and 
P = (— 1) J . In spectroscopic notation, these states are a mixture of 3 (J — l)j, and 3 (J + 1) 7 states. The 
exception is the state with J = 0, for which the orbital angular momentum is a good quantum number. 
Indeed, for J = the function f J ~ 1 (p) does not exist (see Appendix A), thus the function F SlS2 (p) is defined 
only by the second term in (56). Note that ^-mixing states appear for principal quantum number n > 3 only. 



5. The relativistic radial equations for two fermion systems 

It is not possible to write an universal two fermion wave equation, because the adjustable functions have 
different form for different states. Thus it was important to classify all states of the system before deriving 
final radial equations. Now we return to the variational equation (29) and replace the functions F SlS2 (p) by 
the expression (39) for the quasi-states and by (56) for the triplet states. After completing the variational 
procedure we obtain the following results: 

For the states with I = J, P = (— 1) J+1 the radial equations are 

K + n p E) f(p) = ^ [ q2dq JC (p, q) f J (q), (57) 

(27r) J ^UJpUJqilpilq 

where the kernel /C (p, q) is defined by invariant TW-matrices as follows from (29). For the pure singlet states 
with J = 0, according to (42), (48) and (147), the kernel K, (p, q) is 

K,(p,q) = J dp dq(M 12 i2 (p,q) - M 12 2i (p,q) - .M2112 (p,q) + ^2121 (p,q)) , (58) 
For quasi-singlet and quasi-triplet states (J > 1) we have 

JC(p,g)= -i£ f dpdqC { flY 2aia2 M slS2aia2 (p, q )Yr%mj mj ^), (59) 

siS2<7 1 a 2 mj 

where the coefficients Cj S ^^ S2(Tl ' 72 are expressed through the coefficients Cj^™ 3132 

r ,(s)s 1 s 2 vi<J2 _ r ,(s)Jm 3132 (s)Jm„ 1 „ 2 . \ ^ /„(s)Jm n „ 2 \ 2 , . 



v\v 2 mj 



where s = s s , s t . In (60) we have summed over mj, because of the 2J+ 1 energy degeneracy. 

For the triplet states with I = J =F 1, we have two independent radial functions f J ~ 1 (p) and f J+1 (p). 
Thus the variational equation (29) leads to a system of coupled equations for f J ~ 1 (p) and f J+1 (p). It is 
convenient to write them in matrix form, 

m\m2 f q 2 dq 



( LUp + n p -E)F(p) = -^ / 7 == K(p,g)F( g ), (61) 
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where 



and 



F(p) 



f J -\p) 
f J+1 (p) 



K(p,q) 



fcu ip, q) {p, q) 
^21 {p, q) IC 22 {p, q) 

The kernels fCij are similar in form to (59), that is 

Kit (p,q) = I dpdciM SlS2aia2 (p,q)F™— *(p)Y™— (q). 

s 1 s 2 u 1 a 2 m J 

However the coefficients Cj 1 ^ 2 ^ 172 are defined by the expression 

^s 1 s 2 cr 1 a 2 _ ^(ir)<im, 1 , 2 „(fr)i j m„ 1 „ 2 , ^ / ^ l (tr)l i m 3 \ 2 
Jmjij Jmj Jmj I / j \ Jmj J ' 



(62) 
(63) 

(64) 
(65) 



where l\ = J — 1, £2 = J + 1- The system (61) reduces to a single equation for f J+1 {p) when J = 0, since 
f J ~ 1 (p) = in that case. 

To our knowledge, it is not possible to obtain analytic solution of the relativistic radial momentum-space 
equations (57) and (61). Thus one must resort to numerical or other approximation methods. Numerical 
solutions of such equations are discussed, for example, in [11], while a variational approximations have been 
employed in [5]. However, in this paper we shell resort to perturbative approximations, in order to verify 
that our equations agree with known results for positronium to 0(a 4 ). We expect that this must be so, 
given that the interaction kernels (i.e. momentum-space potentials) of our equations involve the "tree-level" 
Feynman diagrams only. 

Beyond 0(a 4 ) our equations are evidently incomplete. One could, of course, augment them by the 
addition of invariant matrix elements corresponding to one-loop Feynman diagrams to the existing Ai- 
matrices in the kernels of our equations. Indeed such an approach has been used in a similar though 
not variational treatment of positronium and muonium by Zhang and Koniuk [14, 15]. These authors 
show that the inclusion of single-loop diagrams yields positronium energy cigenstates which are accurate 
to O (a 5 , a 5 In aj . However such ad-hoc augmentation of the kernels would be contrary to the spirit of the 
present variational treatment, and we shall not pursue it in this work. 



6. The kernels in semi-relativistic expansion and the non-relativistic limit 



For perturbative solutions of our radial equations, it is necessary to work out expansions of the relevant 
expressions to first order beyond the non-relativistic limit. This shall be summarized in the present section. 
We shall do the calculation in the Coulomb gauge, in which the photon propagator has the form [16] 



1 



1 



A>o(k)^, £> 0I (k) = 0, D«(fc")= 



kj 



(66) 



where = (u) p — ui g , p q). 

To expand the amplitudes M of (30), (31) up to the lowest non-trivial order of (p/m) , we take the 
free-particle spinors to be 



u(p,i) 



y 1 + 8m 2 J 

as discussed in Appendix C. Analogously, the photon propagator takes on the form 



(?'P) 

2m 



v(p,i) 



2m 



(67) 



A) (p-q) = - — ^l. £> fe *(p-q)~-- — —2 

(p - q) (p - q) 



Ski - 



(p - q) k (p - q)i 
(p-q) 2 , 



(68) 
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Corresponding calculations give for the orbital part of the .M-matrix 



M ope(orb) / \ 
J SXS2&1&2 \Hi V) 



\2 



J 1 1 / 1 /mi m 2 \ / 1 q-p \ (pxq) 

= ^192 < -2 + o — + 7 + 7 72 + 7 u r <Wi<W2- (69 

i(p-q) m i m 2 \H m 2 rmj yA (p-q) 2 / (p-q) 4 /J 

The linear spin terms, which are responsible for spin-orbit interaction, become 

MZi-l(c l ,p) = -^^ ((^ + 2)^-(^ + 2)^).-^ M , (70) 
4mim 2 \\ m i / \ m 2 / / (p - q) 

Here ~a \ and er \ are spin matrices of first and second particles respectively, which are defined by the 
following operations ~a 1 (p ai x S2 = (V i<£ CTl ) X S2 > ~^2'P ai Xs 2 = L Pa 1 (~^2X S2 )- Thc quadratic spin terms or 
spin-spin interaction terms are 

KA ope(s-s)( x klQl t t J ("?2 • (P- q)) ("CTl • (P- q)) , _> J ,_ 1 , 

Cii«(lP) = i^7<^2 j (p _ q) 2 + ' CT 2)|^ CT1 X S2 . (71) 

The annihilation contribution, which arises for a particle-antiparticle system [10], is given by the term 

A& lCT2 (p,q) = -TT^dxk {^i ' ^2}^X S2 , (72) 

where we have excluded a divergent term, which appears in the Coulomb gauge. Here mi = to 2 = m and 
qi = q 2 = e. 

The kernels are calculated from (58), (59), (64) with (69)-(71). They consist of three parts, namely 

K (p, g) = JC^ (p, q) + /C< s -°) (p, g) + /C^ 8 ) (p, g) . (73) 
The singlet state (^ = J = 0, P = -1, ) 

Details of the calculations can be found in Appendix D. We use the notation: z = (p 2 + g 2 ) /2pg, and 
Q\(z) are the Legendre functions of the second kind [17]. The contributions of the various terms to the 
kernel are as follows: 

Orbital term 

pq 

+ 131S^((^ + ^1 + 1 \ (P + g) Qo{z) + 2Ql{z) _(!!tL + !^ + 2 )) m 

2miTO 2 \\ m 2 mi / \q P J \m 2 m\ J J 

Spin-orbit interaction 

K is9l){s - o) ( P ,q) = 0, (75) 

Spin-spin interaction 

jC {s a i){s-s ){ ) = ^31S1, (76) 

TOlTO 2 

The quasi-states (£ = J(J>1),P= (-1) J+1 , s = 0, l) 

The contributions of thc various terms to the kernel are as follows: 
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Orbital term 



pq 



+ ^ 2 _ mi + rr* _ ( j _ 1} Wp + £ . gj(z) + 2 (J + 1} Qj+i[z) 



2mim 2 \\m 2 mi J \q p 

Spin-orbit interaction 



2mim 2 \ \m 2 mi / \m 2 mi 



(Qj+i (z) - Qj-i (z)) , 



2J+1 



Spin-spin interaction 

KS'^M = C\^- ( P - + q -) Qj (z) -Cl^-^— {JQ J+1 (z) + (J + 1) Qj-i (z)} , 
2mim 2 \q pj mim 2 2 J + 1 

where Ciand C 2 are defined in section 4. 

The triplet states (t = J - 1 ( J > 1), t = J + 1 (J > 0) P = (-1) J ) 

From (64), it follows that the kernels /Ci 2 and /C 2 i are responsible for mixing of states with £ 
and I = J + 1 . These kernels have the form 

K ia (p, q) = (p, g) = " glg2 f ^Q J+ i (*) + (*) - 2Qj (z)) , 

5mim 2 (2J + 1) \5 P / 

where only spin-spin interactions contribute. For kernels ICn and /Cn we have from (64): 
Orbital terms 

Mi (P>«) = — + — +2-J - + - Qj_i(z) 

p<7 2mim 2 \m 2 mi / \9 P/ 

H </<2j \z) - 1 V 2 <5j.i, 

mim 2 2mim 2 \m 2 mi / 

~(o), * 27rqig 2 ^ , s . 7T9i92 /mi ra 2 \ /p g\ 

M2 (P. 9 = <3./+i ( z i + o — + J ) ~ + ~ Qj+iV z i) 

pq 2mim 2 \ TO 2 mi J \q p) 

+ I^i (J + 2)g J+2 ( Z ), 

mim 2 



Spin-orbit interaction 



Ua-o) i \ ""3i92 J - 1 /m 2 mi , 

Mi (P. = o T/ — 7 + + 4 (Qj (z) - Qj_ 2 (z)) 

2mim 2 2J — 1 \mi m 2 



K (s-o) , v 7rgi92 /m 2 mi \ J + 2 
2mim 2 \mi m 2 / 2J + J 



Spin-spin interaction 



2mim 2 2J+l\\g p/ / 
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(P, ?) = o„ i n t i o ((- + -) (*) - 2Q J+ 2 (*) ) . (86) 



J_(Y£ + £ 

2mirn2 2 J + 3 \\g p 
Annihilation term (particle-antiparticle system only, with mi = m 2 = m, q\ = qi = e 



2?re 2 

JC ann { Pl q) = -S JA . (87) 

We note that in the non-relativistic limit the only terms that survive are the first terms of the orbital part 
of the kernels. They have the common form 2Trqiq2Qe(z)/pq, where 



pq f dqdp 
2W (P - qf 



Qt(z) = £ / ^*Y m t (p) Y m ' (q) . (88) 



Thus, equations (57), (61) reduce to the form 



K + n p E) f(p) = Hh^l f dq—f=L=Q e (z)f(q). (89) 



If we expand the relativistic free-particle energy 



1 / P \ \ n I n , 1 / P 



i»i 1 + j r , n p ~m 2 ( l + -(^) ), (90) 



2 \mi J I \ 2 \rri2 

we obtain the two particle Schrodinger radial equation in momentum space [18] 



,2 



T--e)f e (p) = — dqqQ e (z)f(q), (91) 
2m r / ?rp J 

where a = q^l^- The solutions of these equations are well known. They are given in Appendix D (formula 
(201)). 



7. Energy eigenvalues and relativistic corrections to O (a 4 ) for arbitrary mass 
ratio 

The energy eigenvalues E Hi j for the states, when I = J, can be calculated from the equation 

P OO C-OO 

e dpp 2 f J (p)f J ( P ) = / dpp 2 (u p + n p )f J ( P )f J (p) + 

Jo Jo 

(2tt) Jo v^piip Jo 



For the £ = J =F 1 triplet states we have the matrix equation 

' + 



POO POO 

E dpp 2 ¥\p)¥(p) = / dpp 2 (w p + n p )Ft(p)F(p) 
Jo Jo 

-iB r -Pw r ^-Ar^Ft^). (93) 



To obtain results for _E to O (a 4 ) we use the forms of the kernels expanded to 0(p 2 /m 2 ) (expressions 
(74)-(87)) and replace f i {p) by their non-relativistic (Schrodinger) form (201) (see Appendix D). The most 
important integrals, which we used for calculating (92) and (93), are given in Appendix D. In Appendix F 
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we show that the contribution of the kernels /Ci 2 and KL21 in (94), which are responsible for coupling of the 
equations (61), is zero. Thus in the framework of the present approximation, the energy corrections for the 
triplet states with £ = J — 1 and I = J + 1 can be calculated independently. We present the results in the 
form Ae = E + a 2 ^l r — M, M = mi + m 2 , m r = mim 2 / (mi + m 2 ). 

Singlet states (£ = J = 0, P = —1) (which include the ground state) 
The kinetic energy corrections 



Orbital energy corrections 

Spin-orbit energy corrections 
Spin-spin energy corrections 

The total energy corrections 



A (sg )(o) _ (_l + ^L- (Q 5 ) 

£n ~ n 3 \ 2 M Mn)' ( ' 

Ae (sg)(s-o) = 0j ( 96 ) 

Ae ( )(s - s) = _^2_m L ^ 

n 3 M y ' 



n 3 \2 + M 8n I M J ) ' 1 > 



Quasi-singlet, quasi-triplet states {l = J (J > 1), P = {— 1)' 7+1 ) 

The kinetic energy and orbital energy corrections are common for quasi-singlet and quasi-triplet states: 



(fe ) a-m r 1 3m r \ / 2 3 1 

= 2r^~ ( 1 - ~M ) 1 27TT " 4^ ' ' (99) 



Orbital energy corrections 



A£ "> J ~ n 3 M I 2 J + 1 J' ( J 



Spin-orbit energy corrections 



x -C| ((m 1 + m 2 ) 2 +2mim 2 J + 2CiC 2 v/7(jTl) (m| - m 2 

Spin-spin energy corrections 



a (a-o) = ■ _j ■ j 

£ ™< J - 2n 3 M 2 (2J+1)(J+1)J 



A,M _ a4m r m r C l nn .,, 

" n 3 M (2 J + 1) (J + 1) J' 



The coefficients C\ and C 2 are defined here in the same way as in formulas (78), (79). 
Thus, the total energy corrections are 

where upper and lower signs correspond to quasi-singlet s — s s — and quasi-triplet s — s t — I states 
respectively. 

Triplet states (i = J - 1 ( J > 1), P = (-1) J ) 
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Kinetic energy corrections 



Orbital energy corrections 



( t ,)( fe ) = _oSn L ( _J_ __3_\ / 3mA 

n < J n 3 \2J-1 8nJ\ Mj' { ' 



A (tr)(o)_ cAtv / 3^ (5.7,1 l\m r 



Spin-orbit energy corrections 

(tr)(s-o) _ a 4 m r 1-Sjs f 2mA 
A£ "< J - n 3 2J(2J-1) V M J ' UUbj 

Spin-spin energy corrections 

A Jtr)( s - s ) a4fn r ( 1 - 6.7,1 2 \m r 

AS J = n 3 "" U(2J+1)(2J-1) " 3^ J M • (1 ° 7) 
The total energy corrections 

The annihilation term for positronium-like systems is 

A£{ a nh) = e^ 6ji (109) 

Triplet states = J + 1 (J > 0), P = (-1) J ) 
The kinetic energy corrections 

Orbital energy corrections 

/,-)(„, _ n ;;/.. ^ 3 1 ^ //) 

Spin-orbit energy corrections 



Ap - _r_^L _L nil) 

A£ "' J ~ n 3 Uj + 3 n) M 1 [LLL) 



A (tr)(s-o) _ a 4 m r _ 

A£ "' J _ n 3 2(J+l)(2J + 3) V M ' ' 



2m, 



Spin-spin energy corrections 



a (tr)(s-s) = a m r i rrv 

™< J ~~ n 3 (J + 1) (2 J + 3) (2 J + 1) M ' k ' 



The total energy corrections 

. ,„ a A m r ( 1 f I 2m r 1 

A <7 = - 



n" 



3 V 2J + 3 V 2(J+ 1) M 2 J 



tK (>-£))■ <>»> 



For two equal masses our calculations agree with positronium results [19]. In the limit when one of the 
masses becomes infinite, say m 2 — > oo, the above results reduce to those obtained for a one-electron Dirac 
equation in a static Coulomb potential (to O (a 4 )), namely 

a 4 mi ( 1 3 \ ,„ , 

" J = ~ V2JTT ~ 8n) " (115) 
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Indeed, when m 2 — > oo not only the total angular momentum, but also 3\ = (Li + Si) and S 2 are 
independently conserved. Thus, in this case, we can replace the quantum number J by j — 1/2 in (99)- 
(103)(when s = 0) and (110)-(114), and by j + 1/2 in expressions (99)-(103)(when s = 1) and (104)-(108). 
In other words 

J -> j - 1/2 for | sg q ) and \ tr e=J+1 ), (116) 
J -> j + 1/2 for \tr q ) and \tr i=J -i). 

The quantum number j = \£ ± 1/2 1 belongs here to the particle with mass mi. 

The results (98), (103), (108) and (114) agree with the calculations of Connell, based on a quasipotcntial 
reduction of the Bethe-Salpeter equation [20], those of Hersbach, who used a relativistic Lippman-Schwinger 
formulation [21], and those of Duviryak and Darewych based on a two-fermion Breit equation [22]. Grandy 
[8] obtained the same results from a perturbative solution of the coupled, nonlinear Dirac equations (8) and 
(9). 



8. Fine and hyperfine structure. Recoil effects. 

In this section we shall analyze the formulas obtained in the previous section and we shall apply them 
to the energy spectra of some exotic atoms. To compare our calculations with experimental data we shall 
use the standard spectroscopical notation. First of all, it follows from (94)-(96) and (104)-(106) that the 
difference between energy levels of 2S and IS states is given by formula 

E(2S)-E(1S) = ^ + ^ (11 + 15^). (117) 

Note that the formula (117) ignores hyperfine splitting, that is we exclude the spin-spin interaction. The 
fine structure of 2P-state follows similarly from (104)-(106) and (110)-(112), provided that we exclude the 
spin-spin interaction (107) and (113) 

AE fs (2P) ee E (2P 3/2 ) - E (2P 1/2 ) = ^ (l + ^) . (118) 

The HFS of l^i/s and 2S 1/2 is obtained from (98) and (108): 

AE hfs (lS 1/2 )=a i m r ^ f , (119) 



AE hfs (2S 1/2 ) =a 4 m r ^. (120) 

Actually this hyperfine splitting, expressions (119) and (120), arises from the difference of spin-spin terms 
(97) and (107). The formulae (119) and (120) give the usual Fermi splitting [18]. 

The HFS of states with £ > is more complicated. For each I > we have two sorts of states with 
different J and s: states with J = £+ 1 and J = £, s = s s and states with J = £ — 1 and J = £, s = s t . 
From (103), (108) and (114) we obtain the general HFS formulae for all quantum numbers n, I and for any 
mass ratio. 

* j~i ( fi \ a tv a ^m r 1 /2^+1-r 1 2m r 1 \ 

AE hfs (n,£,s s ) ee A< J= , +1 - Ae n ,j^ 3 = [ + -Jf^) , (121) 

4P , . , . . f r a 4 m r 1 /21+1-C 1 2m r 1 \ /ioo ^ 

AE hfs (n,£, st ) ee Ae nM - Ae* M = { u{£+1) + ^^l) ■ ^ 

The quantity £ is defined by (53), but with the quantum number J replaced by £. 
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For the particular case when m 2 >> mi we obtain from (121) and (122) 

. , „ a A mi 8 (£ + 1) mi , „„. 

AE hfs (n,£,s s ) = — ^- i-= — '- -, 123 

. ^ , a 4 mi 8£ mi , 

AE hfs (n,£,s t ) = — ^- o -. 124 

1 V ; (2^+l) 2 (2£- l)m 2 ' v ' 

to leading order in mi/m 2 . 

The approximate results (123) and (124) are the HFS formulae usually quoted in the literature (e.g. [18] 
p. 110, and [24] p. 836). As mentioned in section 7, in the one-body limit (777,2 — > 00), states with labels 
s s and s t can be characterized by the quantum numbers j = £ + 1/2 and j = £ — 1/2 respectively. In this 
one-body picture, the labels s s and s* should be replaced by j in the formulas (121), (124). Note here that, 
because of the mixed nature of the states n£i-i/ 2 {£ — J, s t — 1) and n£i + i/ 2 {£ = J, s s = 0), the spin-orbit 
interaction contributes to the spliting as well. The splitting is largest when mi — m 2 and disappears for all 
states in the one-body limit. 

It is easy to get the HFS of 2P 1/2 , 2P 3/2 states from (121) and (122): 

Aa,.(2P V3 )= ( ,^(i-g + i^), ,125, 
A B „.(2P 3/2 )=«V(^-g + 4^), (126, 



where f -1 = ^8 (^f^) 2 ~. Note that, if we use the pure triplet state instead of the mixed states, 
we get the following expressions for hyperfine splitting: Ehf s (2Pi/ 2 ) = a 4 m r + jttjj), Ehf s = 
a 4 m r (<^g + ■ These give the correct results only for mi = m 2 (positronium). 

The results of our calculations and comparison with experimental data for muonium, hydrogen and 
muonic hydrogen are summarized in Tables 2-4. In each table the first two entries reflect, primary, fine- 
structure splittings, while AE^ corresponds to the O (a 4 ) recoil contribution to the Lamb shift, E (25 ly / 2 ) — 
E(2Pi/ 2 ). Note that except for the hyperfine splitting of the ground state, which has been measured very 
accurately, the experimental results are quoted in the tables to the number of significant digits available. 
The calculated results, we must remember, are accurate to O (a 4 ) only, where as the experimental results 
are exact except for experimental error, that is they contain all orders in a. 

Although hydrogen and muonic hydrogen are not strictly leptonic atoms, we can take into account the 
proton anomalous magnetic moment in the hyperfine splitting in the same phenomcnological manner as was 
done in [23]. Thus, to get the corrected result for HFS in hydrogen and muonic hydrogen, we multiply the 
expression for AEhf s by the factor (1 + fc p ), where k p = 1.792847 is the anomalous magnetic moment of 
the proton. The last columns in tables 3, 4 give the HFS corrected by the factor (1 + k p ). The anomalous 
magnetic moment of the muon is very small, nevertheless for the sake of completeness, we provide the same 
corrections for muonium, where we use the factor (1 + a^) (a M = 0.001166) to get the result of the last 
column in Table 2 from those of the third column. We note that doing so improves the agrement between 
the calculated and the experimental HFS results in most cases, and dramatically so for hydrogen and muonic 
hydrogen. 

Table 2. Energy shifts, fine and hyperfine structure in muonium for 77 = 1, 2. 
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e 




experiment [24] 


theory 


a^-correction 


E{2S 1/2 ) - 


- E(lS 1/2 ) 


2455. 3 1 Hz 


2455.54 1 Hz 




AE fs (2P) 


10922 MiJz 


11000.85 MHz 




AE L 


1047 MHz 


69.59 MHz 




AEhf s 


^1/2 




4463 Mi? z 


4453.84 MiJz 


4459.03 MHz 


AEhf s 


^1/2; 




558 MHz 


556.73 MHz 


557.38 MHz 


AEhfs 


2^1/2, 




186 Miiz 


185.80 MiJz 


186.01 Miiz 


AEhfs 


2^3/2, 




74 MiJz 


74.43 MHz 


74.52 MiJz 



Table 3. Energy shifts, fine and hypcrfine structure in hydrogen for n = 1, 2. 



H 


experiment [25], [26] 


theory 


/cp-corrcction 


E (2S 1/2 ) - 


- S(15i/ 2 ) 


2466 THz [25] 


2469.99 THz 




A£ /s (2P) 


10969.04 Mi? z [25] 


10954.77 MHz 




A£ L 


1057.85 MHz [25] 


7.94 MHz 










1420.4058 MHz [26] 


508.0265 MHz 


1418.8403 Miiz 




v 2Si/2 




177.5569 Miiz [25] 


63.5033 MHz 


177.3550 MHz 




2^1/2, 




59.1721 MHz [25] 


21.1698 MHz 


59.1241 Mifz 




2P3/2 




23.6541 Miiz [25] 


8.4696 MHz 


23.6543 Miiz 



Table 4. Energy shifts, fine and hypcrfine structure in muonic hydrogen for n = 1, 2. 



M -ii 


theory [27], experiment [28] 


theory 


fcp-correction 


£ (25 1/2 ) - E(lS 1/2 ) 




1.90 




AEfs (2P) 


8.35 meV [27] 


9.95 meV 




A£ L 


-202.07 meF [27] 


+ 1.02 meF 




A£7 fc/8 (lS 1/2 ) 


183 meV [28] 


65.33 meV 


142.36 meV 


A£7 fc/8 (2S 1/2 ) 


22.75 meF [27] 


8.17 meV 


22.81 meV 


A^ /8 (2P 1/2 ) 


7.96 meV [27] 


2.79 meF 


7.79 meV 


A^ /8 (2P 3/2 ) 


3.39 meV [27] 


1.16 meV 


3.23 meV 



From tables 2-4 we see that the O (a 4 ) calculated results agree with observation remarkably well, except, 
of course, for the Lamb shift, which is dominated by corrections of higher-order that a 4 . The present theoret- 
ical results for AEl reflect only the O (a 4 ) recoil effects, and these provide a small, though not insignificant, 
contribution. Note that if we use the approximate mi/m 2 expansion formulas (123) and (124), for example for 
HFS of 2P1/2 and 2P 3 / 2 states we get the result (with a M and k p corrections) AEhfs (2P1/2) = 188.50 MHz, 
AEhfs (2P3/2) = 75.40 MHz for muonium, and AE hfs (2P1/2) = 10.47 meV ', AE hfs (2P 3/2 ) = 4.19 meV 
for muonic hydrogen. These approximate results do not agree as well with observation as the values calcu- 
lated without the m\jm 2 expansion, and hence they differ from what we give in last columns of tables 2 and 
4. Note that this difference can be of the same order as the contribution of higher order in a for HFS. The 
mi/m.2 expansion approximation is not significant for hydrogen for which the mass ratio is very small, but 
it is appreciable for muonium and particularly for muonic hydrogen. 

Since our results arc true for arbitrary mass ratios we can speak here about recoil effects rather than 
about recoil corrections. The recoil effect contributes to the Lamb shift (shift between 2,Si/ 2 and 2Pi/ 2 
levels). However the states 2S' 1 / 2 and 2P 1 / 2 are split (hypcrfine splitting as discussed below) because of the 
spin-spin interaction, so to calculate the Lamb shift we have to take the energy corrections for each state 
without spin-spin terms. For the state 2S\/ 2 it is the sum of the energy corrections (94) and (95), while for 
the 2Pl/ 2 state we take the sum of (110), (111), and (112). Thus, the formula for the O (a 4 ) contribution 
to the Lamb shift produced by the recoil effect is 

AEt l) = E (2S 1/a ) - E (2P 1/2 ) = (127) 
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It is easy to see that the shift (127) is largest when m\ = m 2 , and disappears in the one-body limit, when 
m 2 — > oo. 

There is a well known formula for the recoil correction of order a 4 obtained on the basis of the Breit 
interaction [29], [30] 

^-^(jT^-Trb)' 1 -'")- <128 » 

which predicts the following contribution to the Lamb shift: 

. P (rc() A a a A m r to 2 

AE y L = Ae 2 ,o,i/2 - A£r 2 ,i,i/2 = 73 £• 129) 

' ' 48 m 2 

Comparison of the formulas (127) and (129) shows that the total recoil effect (127) is larger in magnitude than 
(129), and has opposite sign. The difference between (127) and (129) stems from the spin-orbit interaction: 

If we ignore the term + 2^ er 2 in expression (70), we get instead of (106) and (112) the following result 
for the spin-orbit contribution: 
for i = J - 1 

A (tr)(s-o) _ a 4 m r 1 - S JA f m 2 r \ 
A£n > J ~ n3 2J(2J-1)\ 1 ^J' 

for £ = J + 1 

A (tr)(s-o)_ a 4 m r 1 / m 2 r \ 

A£ "' J _ n3 2(J+l)(2J + 3) y- ml)- 

The relevant terms are the second ones in (130), (131). We can replace the "J"-quantum numbers in (130) 
and (131) using (116): 



1 2 \ / 1 1 



J(2J-1) \J 2J-1J + £+1/2 

and 



(132) 



1 = J: L_ _ _J L_ (133 ) 

(J+l)(2J + 3) J+l 2J + 3 j + 1/2 £+1/2' v y 

It is now easy to see that the second terms in (130) and (131) give the formula (128) which, therefore, 
must be regarded as incomplete. As already mentioned, the proper recoil contribution to the Lamb shift of 
muonium, hydrogen and muonic hydrogen, based on our result (127), is given in Tables 2-4. 



9. Concluding remarks 

We have used the variational method within the Hamiltonian formalism of QED to derive relativistic 
momentum-space wave equations for two-fermion systems like muonium. The trial states are chosen to be 
eigenstates of the total angular momentum operators J 2 , J3 and parity, as well as charge conjugation for 
particlc-antiparticle systems. A general relativistic reduction of the wave equations to radial form is given for 
arbitrary masses of the two fermions. For given J there is a single radial equation for total spin zero singlet 
states, but for spin triplet states there are, in general coupled equations. We have shown how classification 
of the states follows naturally from the system of eigenvalue equations (36), given our trial state. 

It is not possible, as far as we know, to obtain analytic solutions of our relativistic radial equations nor 
the resulting eigenvalues of the two fermion system that they describe. However, it is possible to obtain 
0(a 4 ) corrections to the energy eigenvalues analytically for all states using perturbation theory. 

We have compared our calculated results with experiment for fine and hyperfine splitting of low-lying 
levels in muonium, hydrogen and muonic hydrogen. We find good agreement for muonium, as well as for 
hydrogen and muonic hydrogen provided that we take into account the anomalous magnetic moment of the 
proton in the later two cases. 
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The method presented here can be generalized to include effects higher order in alpha by using dressed 
propagators in place of the bare propagators. This shall be the subject of a forthcoming work. 
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Appendix A: Total angular momentum operator in relativistic form 

The total angular momentum operator is denned by the expression 

J = y rf 3 xV f (x) (L + S)^(x) + J d 3 x0 f (x) (L + S)<f> (x) , 



(134) 



where L is the orbital angular momentum and S - the spin operator: L = xxp and S = \ ~~o . We use the 
standard representation for the Pauli matrices 



a 
~a 



' 


1 " 




' -i~ 




' 1 


1 





cr 2 = 


i 




-1 



(135) 
(136) 



Using the field operators tp (x) and (f>(x) in the form (17), (18), after tedious calculations we obtain the 
expression for operator J. It consists of three parts: two total angular momentum operators of each particle- 
antiparticle system [10] and the following part, relevant to our considerations 



(137) 




Here L 9 is the orbital angular momentum operator in momentum representation: 

(L q )i = L ql = -i (qxV 9 ). . 



(138) 



Note that these expressions are valid for any t, since the time-dependent phase factors of the form e lulqt 
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cancel out. For the operator J 2 = j\ + j| + j| we have 



where 



+ 



d 3 q 



2L„ 



+*U b *iKi D *+ b *M D ii D «i 



+ ± (bUb a ^DlD 



+L 



^M D l<i D ^ + b *M D l'i D <*'i 



q'+ 



+ (W 3 + L q3 ) (b^b^D^D^ - bl^^D^ 
\ - fata - L q3 ) fe T Vt^Au - b^b^D^A ) 



Lq+ — L q l + iL q2 , 



Lq- — L q l — lL 



(139) 



(140) 



The formulae (137) and (139) apply to the particle-antiparticle system [10] as well if the operators and 
D are formally replaced by S and d. 

The requirements (36) with trial state in the form of (24) lead to the system of equations (F SlS2 (p) = 



IF, 



L3F12 
L3F21 

(^3 — F22 



mjFn 

m,jF 12 
mjF 2 i 

mjF 2 2 



(141) 



( J( J + 1) - L 2 - 2 - 2L 3 ) F n = L_{F 12 +F 21 ) 

(j(J+l)-L 2 -l)F 12 = F 21 +L + F U +L_F 22 
(j(J+l)-L 2 -l)F 21 = F 12 +L + F 1X +L-F 22 

( J( J + 1) - L 2 - 2 + 2L 3 ) ^22 = £+(^12+^21) 
After substitution of the functions F SlS2 , Eq. (38), into the system (141) and (142) we get 

m 12 = m 2 \ = mj, ran = raj - 1, m 22 =mj + l, 

in = £22 = £12 = £21 = e, 



(142) 



(143) 
(144) 
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and 



(J(J + l)-£(£+l)-2m J ).f 11 (p) = ^(£-mj + l)(£ + mj)f( 2 (p) 

W(t~™j + m + mj)fiM (145) 

(j(j + i)-e(£ + i)-i)f( 2 (p) = f 21 { P ) 

+ y/{l + mj){l-mj + \)fii ( P ) 
+ ^(£-mj)(e + m J + l)f 2 ! 2 (p) (146) 
{J{J+l)-i{l+l)-l)f 21 {p) = f 12 (p) 

+ ^(£ + mj)(£- mil + l)fn (p) 
+ V(e-mj)(£ + mj + l)f^ 2 (p) (147) 
(J(J + l)-£{£+l) + 2mj)f 22 {p) = ^/(£ + mj + l)(£-mj)f 12 (p) 

+ V(t + rnj + l)(£-mj)f' 1 (p) (148) 

The solution of this system leads to two categories of relations among the functions /f lS2 (p). The first 
category, which we call the trivial one, is obtained when fd (p) = f 22 (p) = 0. In this case, as it easy to 
see from (145) and (148), we get /f 2 (p) = —f 21 (p). This solution corresponds to the singlet states of the 
system with £ = J (J > 0) (as follows from (146) or (147)). This simple relation allows us to write the 
general formula for the components f% l82 (p) in the following form 

f SlS2 (p) = C^ m ^ f J (p), (149) 

where the radial function f J (p) is common for all components, and the coefficients C^ 9 ^" 1 ^^ have the 
following properties 

(j(sg)m ll _ (j{sg)m 2 2 _ Q (j(sg)m 2 i _ _(j(sg)m 12 _ ^ (150) 

The formulas (149) and (150) with the notation (p) = fstsY (p) make cq. (41) evident. 

The second category corresponds to the triplet states. First of all, it is not difficult to see from (146) and 
(147) that f( 2 (p) = f 21 (p) = f i (p), and after some simple calculations we get: 

for I = J - 1 (J > 1) 



(J-mj)/,'- 1 ^) = ^(J-mj)(J + m J -l)f J - 1 (p) (151) 
(J + mj)/^- 1 ^) = y/(J + mj) (J - mj - l)f J ~ 1 (p), (152) 



for £ = J (J > 1) 



mjfA(p) = -V(J + mj) (J-mj + l)f J (p) (153) 
m.j.fi 2 {p) = y/(J-mj) (J + mj + l)f J (p), (154) 



for £ = J + 1 (J > 0) 



(J+l + mj)f^ +1 (p) = -y/{J - mj + 2) (J + mj + l)f J+1 (p) (155) 
(J+l-mj)f 2 T 2 +1 (p) = -v/(J - mj + 1) (J + mj + 2)f J+1 (p). (156) 

It is convenient to introduce the table of coefficients Cj^ ms , which represent the relations in (151)-(156). 

Table 5. The C-G coefficients for triplet states (total spin 5=1). 
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J - 1 



J 



l = J + l 



m s = + 1 

(J+m.,-l)(J+mj) 
.7(2.7-1) 



TJT 



(J+mj)(J-mj + l) 
7(7+1) 



TJT 



(,J-mj + l)(J-m,j+2) 
(.7+1) (2. 7+3) 



TJT 



m s = 

' (J-rnj)(J+rnj) 
J(2,7-l) 



T72- 



(./(■/+D) 1/2 



/ (J-m,,+l)(J+m,,+l) \ 
1, (,7+l)(2J+3) ) 



TJT 



(J-mj-l)(J-mj) 
,7(2 J— 1) 



T72- 



(.7-m. f )(.7+mj + l) 
./(./+!) 



T72- 



/ (J+m,,+2)(J+mj+l) \ 
1, (,7+l)(2.7+3) ) 



TJT 



Thus, we can write the relations between the components fl a (jj) in the compact form 

Cp) = 45? b -/ < (p)- 



The coefficients Cj^^ ms coincide with the usual Clebsch-Gordan coefficients for total spin 5=1 except for 
a factor 2 in the denominator. The expression (38) can now be written in an explicit form 



(157) 



f s1S2 (p) - c^/- l ^f-\ P )Y^r (p) + cf2Y )ms f J tp)Y? slS2 (p) + C +1)m */ J+1 (p)C t2 (p)- 

(158) 

However, as is shown in Appendix B, the first and third terms have parity P = (— l)' 7 , while the second 
term has parity P = (— 1)' /+1 . Thus, we get the result (56) by suppressing the second term in (158). The 
second term in (158) is associated with the singlet solution (42) for the mixed-spin states, which have the 
same parity. 



Appendix B: Parity and charge conjugation 

We consider the application of the parity operator to the trial state (24): 

V Atrial) - E / d 3 P F SlS2 (p)Pbl Si Dl pS2 |0> = W d^Fs^Vbl^VDl^V-'V |0> . (159) 

sis 2 sis 2 

Making use of the properties 

Hs^' 1 = V P bl psl , VDl^V- 1 = -V P Dl 2 , V |0) = |0) , (160) 
where i] p is the intrinsic parity ((f] P ) 2 = 1), it follows that 

V \Arial) = E/ d 3 P F SlS2 (p)VblDl pS2 |0) = -E / d 3 pF SlS2 (-p)bl s Dl pS2 |0) 

sis 2 sis 2 

= PE / rf 3 P F SlS2 (p)6t sijDl . pS2 |0), (161) 

sis 2 J 

where the parity eigenvalue P depends on the symmetry of F SlS2 (p) in different states. For the singlet states 
(£ = J) we get from (36) F SlS2 (-p) = (-1) J F SlS2 (p), so P = (-1) J+1 . For the triplet states with £ = J 
we get from (47) F SlS2 (-p) = (-1) 3 F SlS2 (p), so P = (-1) J+1 . For the triplet states with £ = J± 1 we get 
from (56) F SlS2 (-p) = (-1) J+1 F SlS2 (p), so P = (-1) J . 

Charge conjugation is associated with the interchange of the particle and antiparticle. Applying the 
charge conjugation operator to the trial state of the form (24) we get 

C |VW*> = E/ d 3 pF SlS2 (p)Cbl Si dl pS2 |0> = E / d^pFs^CblJ-'Cdl^C-'C |0) . (162) 

sis 2 sis 2 
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Using the relations 



Cb^C- 1 = V C 4 Sll Cdl p J-' = r, c tf 



C |0}= |0). 



where (i] C ) 2 = 1, we obtain 

C Atrial) = E / d 3 pF SlS2 (p)Cbt Sl dl pS2 |0> = -E/ d 3 pF S2Sl (p)bl Sl dl pS2 |0) 

SiS 2 S!S 2 

= c£/ rf 3 pF SlS2 (p)fet si ^ pS2 |0), 



(163) 



(164) 



where the charge conjugation quantum number C depends on the symmetry of F SlS2 (p) in different states. 
For the singlet states (I = J) we get from (42) F SlS2 (-p) = (-1)' /+1 ,F SlS2 (p), so C = (-1)' 7 . 
For the triplet states with £ = J we get from (47) F S1S2 (— p) = (— 1) J F S1S2 (p), so C = (— 1) J+1 . 
For the triplet states with £ = J ± 1 we get from (56) F S1S2 (— p) = (— 1)' 7+1 F slS2 (p), so C = (— l)' 7 . 



Appendix C: Expansion of the spinors and .M-matrix elements 

We recall the form of the particle spinors: 



u(p,i) = N p 



1 



where 



¥>2 = 



lo p + mi 



2toi 



The antiparticle or "positron" representation for the Vi(p) spinors has the form 



v(p,i) = N p 



Up+mi 
1 



X.i 



where 



The normalization is 





' " 




' 1 " 


Xi = 


1 


X2 = - 






u(p,i)u(p,j) = Sij, v(p,i)v(p,j) = -Sij. 
Expanding in powers of p/mi and keeping the lowest order non-trivial terms, 

(~a-p) (~a-p) 



cj p + mi 2mi 



we obtain the result 



v(p,i) ~ 1 



8m? 



w p + mi 




2mi 


1 




(V-p) 


<^ = 


2mi 




r (^-p) ~ 




2mi 

i 


Xt = 



\ T 8m- 

(^■P) " 

2mi 



2mi 



Pi. 



(165) 
(166) 

(167) 

(168) 
(169) 

(170) 
(171) 

(172) 
(173) 
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The A4- matrix elements (69)-(72) have the following nonzero components: 
Orbital 



Spin-orbit 



Spin-spin 



= MmW = M222W 



= m<i2 



1 , 1 ( m 2 | mi A _ ((p-q)'P) 
p-q) 2 8mirri2 l^mi m 2 J 

( U12. _L m l ^ 



mim 2 (p-q) 



mim 2 (p-q) 

qp 

2m 1 m 2 (p— q) 2 



m 2 mi 
— + — +4 



(s-o) 



(q»p) = 



8mim 2 V mi m 2 



+ 4 



(p-q) 



Mn 2 o) (q,p) = M s i22 o; (q,p) = 



(s-o). 



«<7l<72 



(s-o) 



(q,p) = A^ 2 2 21 o) (q,p 



4mim2 \ m 2 
igig 2 / mi 



(p-q) 



•Wn2i o) (q.p) = xi 2 22 0) (q:P) 



4mirri2 \m 2 
%<?2 



+ 2 



Mm 5 (q. p) = ^2212 J (q> p) 



(s-o). 



4mirri2 \ TO i 
iqiq 2 ( m 2 



^ + 2 



^i S 272 o) (q,p) = -Mi2r(q,p) = 



(s-o)/ 



4TOiTO2 Wl 



«<7i<72 / m 2 



+ 2 



8miiri2 \mi mi 



(p- 


q) 






(p- 


^2 ' 

q) 




P-93 


(p- 


q) 


P+<?3 - 


P3<7+ 


(P 


^2 ' 

q) 


\ P-<7+ 


-p+q- 


/ (P 


-q) 



MmW) - M222W) - -^m 2 s) (q,p) = (MmW))' 

*9l<?2 / (?3-P3) 1 \ 



4mim 2 \ (q - p) z 



MmW = (4i1q,P)) , = -^;n ) (q,P))' = -4m ) (q 1 P 
MmW) = (^^2 ) (q,P))* = -^&2 ) (q,P) = -(-M^m ) (q,P)) 



4mim2 



(q-py 



M*i22 a) (q,p) = -(M a 2 n ) (q,p))* = 



»gi<? 2 g_-2g_p_+pg 
4mim 2 (q-p) 2 



.,(s-s)/ \ a ,,(«-*)/■ \ ^192 I q+q- - p+q- - p-q+ + p+p- 

■^1221 (q.p) = -^2112 (q»p) = I =5 2 



4mim 2 



(q-py 
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Here p = (p 1 ,p2,p 3 ) and 

P+ = Pi + iP2 = {0, ip) ■ 

P- = Pi- iP2 = \J^-pY{- 1 (0, if) , 

P3 = ^— P Y { (9^). 
For a particle-antiparticle system the annihilation components are 



kAanh kAanh ^ 

- vt llll — • /vl 2222 — 2ro2 ' 



2 

kAanh kAanh kAanh kAanh ^ e 

Vl 1212 — ■ /Vl 1221 — • /Vl 2112 — ■ /W 2121 ~ ^2 ' 



Appendix D: Some useful expressions, identities and integrals 

The following expressions and identities arc useful for evaluating the M. -matrix: 



1 2^ 



(q-p) 2 IpIIq, , 

A m\ = — \ 

where z = (p 2 + q 2 ) /2pq, and Qa(^) is the Legendre function of the second kind of order A [17]. 

((p - q) • p) 2 = P 2 _ (P x q) 2 
(p - q) 4 (p - q) 2 (p - q) 4 

The angular integration in (58), (59), (64) involves the following integrals 

J dp dq F (p • q) Yy 1 {WT J * (p) - 2tt W ro>l7 (p • q) F (p • q) F/ (p • q) , 

f (pxqf n (J + 1HJ + 2) J(J-l) n , v 

i (P ' ^ Jp~*f (P ' ^ = 2(2J+1) <*> - 2(2J+1) Q (Z) ' 

Here F (p • q) is an arbitrary function of p • q, Pj (x) is the Legendre polynomial of order J. 
The integrals in the form 

J dpc«(p)F;hp)y;"(p) 

can be calculated using the Wigner-Eckart theorem [17]. 

The calculation of the relativistic energy corrections involves the integrals 



o Jo 



oo />oo 



JO 



dp dq pVf J (p)f J (q) = 2n(^-) 5j, , 
dpdqpqf J (p)f J (q)Q J (z 1 ) 
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/ / d P dqpWf J {p)f J {q)Qj{ Zl ) = 

Jo Jo 

/ / dpdqp 3 qf J (p)f J (q)Qj(z 1 ) = tt (^L 
Jo Jo v n 



(198) 



o Jo 



2J+1 n, 



OO /'OO 

2 2 fJ(^\ fJ 





'2 


3 









jT jT W.VW<«»,« M = • (^) 3 - i) ■ (20») 

Here f J is the nonrelativistic hydrogen-like radial wave function in momentum space [18] 

f j ( P ) = f(P) r 2 (w - j - i)! V /2 " j+ v 22(j+i) y + i (2 oi) 

tW-toW U (n+J)! J (nV + l) J+2 ^-^UV + iy' ( J 

where Gn-~j-\ ( x ) are the Gegenbauer functions. 

Appendix E: Diagonalization of the expectation value of the Hamiltonian 

The matrix representation of the perturbing Hamiltonian AH — H — Hnr — M, in the basis of the states 
\sg},\ tr) (48), (49), is (J + 0) 

(V| AH |V) t=0 = [ aij ] {i, j = l,2), (202) 



where the matrix elements 



an - <„| AS | S9 > _ (^ _ ( 3 _ _L) , (203) 



/ i a ft i > a 4 m r ( 1 1 / m r \ 1 \ ,„„,s 

„ 23 = W A ff W = (— - 2J(J + 1)(2J+1) " (3 - ^) , (204) 

« - a 21 - AS |«r) . « r | AS N ) . A^O^_L_. (205) 

Note that, in the case of positronium, the elements an and 022 give the energy corrections for pure singlet 
and triplet states respectively. Diagonalization of this matrix leads to (103) with eigenvectors (54). 

Appendix F: /C12, /C21 kernels for l=j^l states 

The contribution to the energy correction due to the kernel /C12 is 

J dpdqp 2 q 2 IC 1 2( P ,q)f J - 1 (p)f J+1 (q), (206) 

where 

lC l2 {p,q)= CjZ™ f d^d^MZs^A^) Y T;rmT-r\v). (207) 



(JiO- 2 SiS2 
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This requires the following integral 

E CjZi? 2 I d 3 pd\f J -\p)Y^r\p)M7j£-°l (p,q)/ J+1 (?)F™r 2 (q)- (208) 

a 1 a 2 s 1 s 2 

We calculate this form in coordinate space. The Fourier transforming of M SlS2!7ia2 (p, q) is 

M SlS2ai(T2 (p,q) = J d 3 rd 3 r'M SlS2aia2 (r, r') e -*(P-i)-('-'') , (209) 

where the .M SlS20 - 1(72 (r, r') matrix in general is a local operator [16] 

M SlS2<Jl a 2 (r, r') = .M SlS2<TlCT2 (r) S (r - r') . (210) 

We apply this transformation to the Msitia^al (p, q)-matrix (see (71)). Because of the angular integration 
in (207), only the first term of (71) contributes. The Fourier transformation of that term is 

("ct 2 • (p-q))("ai • (p-q)) ("a 2 -r) (gi -r) 

1 o R ' V ' 



4 (p-q) 2 167rr 5 



Furthermore, 



y d 3 P / j - i (p)r J m 7-*(p) e -p- r = i?^v)^rr 2 >), (212) 
| rf 3 q/ j+i ( 9 )r; n ;r 2 *(q) e - jq - r = R j n + \r)Y^r w, (213) 



where 



n / 1 \n 



The associated Laguerre functions (p) are related to the confluent hypcrgcomctric functions by 

L " (P) = { ~ ir »l(\-ri\ F (_A + M "" + 15 P) ' (215) 
The generating function for the Laguerre functions is 



hence 



E <w2* a / ^P^q/ 7 - 1 ^)^-! 182 *^)^^;^ ( P ,q)/ J+1 (5)^r 2 (q) 

0"l<T2SlS2 



aia 2 sis 2 

3a 
167rTOim 2 

x 



J drr 2 ±R J n -'(r)R J n + 1 (r) 



E ^-r/^^ ro -\ iS2 *(?)(^2-?)(^i-?)0 S2 ( ? )- 



<Ti<T 2 SlS2 
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We can show, that 

1 y/JJJ+T) 



but 



crio- 2 sis 2 

/■OO 1 

/ drr 2 -^- 1 (r)i?;{+ 1 (r)=0. 
Jo ' 



'0 

The last expression can be proved in the following way. Let us consider the more general case 

p 00 

/ drrP +2 Ri(r)RUr). 
Jo 

The generating function for R l n (r) is 



2 / (n-l-1)! _ r/n (2rV n+1 u 2 ^ 1 / u 2r 

G ne (r,u) = --^J— — — — T e ' — (-1) ^r^exp^ 



n 2 y ((n + £)!) 3 VW (1 - u) 2£+2 I 1-un 

Then we consider an expression 



/>OC 

/ drrP +2 G nt {r,u)G nt ,(r,v) 
Jo 



X 



It is well known that 



Hence 



[°° drr W± l (n-i-l)l(n-e^ c _ 2r/n (2r\^' 
Jo « 4 V ((n + £y.f((n + l')\f \nj 

u 2t+i v 2f+i f / U v \ 2r\ 

X (1 - uf +2 (1 - vf' +2 ° XP I ll^ + T^J ~) 

A_ \ {n-l-l)\{n-l' -\)\ u 2i+i v 2t+i 

((n + ey.f((n + £>y.f {i- u ) 2£+2 (i-v) 2e ' +2 x 

/ dpp /3 e-'' = r (0 + 1) . 
Jo 



00 , y+2+t+t , , x 2r 

dr — exp - 1 + + — 

\ n J I \ 1 — u 1 — v J n 



and 



/ rfr^ +2 G^(r,n)G^(r,t;) 
Jo 



2- /3 - 1 (n-l-l)\{n-£' 

x 



((n + £)!) 3 ((n + ^')0 3 



(1 — uv) 



J+3+I+I' Tifl + Z + t + t). 
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We expand this expression in a series: 

drrP +2 G ne (r, u) G nt > (r, v) = ^ C w (n, [3, £, £') vPvP' . (226) 

It is not difficult to show [31], that the coefficient C n +i tn +i' are given by the integral 

/>oo 

C n+t , n+t . (n, [3, £,£')= / dr/ )+2 R e n (r)R^(r). (227) 
Jo 

Simple but tedious calculations show that this coefficient is zero for (3 = —3, £ = J — 1, £' = J + 1. Thus 
the kernel IC12 does not contribute to the energy corrections at O (a 4 ) . The same result is obtained for the 
kernel A^2i- 



L 
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